HW 8_Onyango
Brenda Onyango
10/12/2021
Load needed packages
library(tidyverse)## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.1 ──## ✓ ggplot2 3.3.5 ✓ purrr 0.3.4
## ✓ tibble 3.1.5 ✓ dplyr 1.0.7
## ✓ tidyr 1.1.4 ✓ stringr 1.4.0
## ✓ readr 2.0.2 ✓ forcats 0.5.1## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## x dplyr::filter() masks stats::filter()
## x dplyr::lag() masks stats::lag()library(janitor)##
## Attaching package: 'janitor'## The following objects are masked from 'package:stats':
##
## chisq.test, fisher.testlibrary(rstan)## Loading required package: StanHeaders## rstan (Version 2.21.3, GitRev: 2e1f913d3ca3)## For execution on a local, multicore CPU with excess RAM we recommend calling
## options(mc.cores = parallel::detectCores()).
## To avoid recompilation of unchanged Stan programs, we recommend calling
## rstan_options(auto_write = TRUE)##
## Attaching package: 'rstan'## The following object is masked from 'package:tidyr':
##
## extractlibrary(bayesplot)## This is bayesplot version 1.8.1## - Online documentation and vignettes at mc-stan.org/bayesplot## - bayesplot theme set to bayesplot::theme_default()## * Does _not_ affect other ggplot2 plots## * See ?bayesplot_theme_set for details on theme settingExercise 1
- The steps for grid approximation are:
- Define a discrete grid of possible theta values
- evaluate the prior pdf of f(theta) and liklihood function L(theta|y) at each theta grid value 3)Calculate the product of f(theta) * L(theta|y) at each theta grid value and then normalize the products so they sum to 1
- randomly sample N theta grid values
- To make the approximation more accurate I would change the first step when defining a discrete grid. Increasing the number of grid values creates approximations that better match true posteriors.
Exercise 2
#Using code Yue shared to include picture
knitr::include_graphics("/Users/brendaonyango/Desktop/IMG_4318 copy.png")
Exercise 3
If the chain is mixing too slowly it will have only explored a limited range in the first several thousands iterations. It will overestimate the plausibility of pi values in this range and underestimate the plasuibility of values outside the range.
The chain having a high correlation will not estimate the correct peak, if there is one, in the posterior.
Chains that get stuck are overestimating in the values where they are stuck and will produce peaks in the posterior that are not there.
Exercise 6.4
MCMC diagnostics are important because simulations are not perfect and diagnostics, used holistically, can identify how to improve a MCMC chain.
MCMC simulations are helpful because they are an efficient and flexible alternative to grid approximations. Grid approximations for more complex problems take a long time and grid approximations are not good for estimating posteriors when there are multiple dimensions. MCMCs can.
Rstan combines R with the Stan “engine”; Stan is written in C++ and is good for Bayesian modeling.
How to put Markov chains in applied language for social problems/questions.
Exercise 6.5
We’re given Y|pi~Bin(n, pi) and pi~Beta(3, 8) and n = 10 with Y = 2.
- Use grid approximation with values pi E{0, 0.25, 0.50, .75,1} to approximate posterior
#step 1: defining a grid with 5 pi values
grid_data <- data.frame(pi_grid = seq(from = 0, to = 1, length = 5))
# step 2: evaluating prior and liklihood at each pi
grid_data <- grid_data |>
mutate(prior = dbeta(pi_grid, 3,8),
liklihood =dbinom(2,10, pi_grid))Next is approximating the posterior using the product of the liklihood and prior and then normalizing them.
#step 3: approximate the posterior
grid_data <- grid_data |>
mutate(unnormalized = liklihood * prior,
posterior = unnormalized/sum(unnormalized))
#confirming that posterior approximation sums to 1
grid_data |> summarize(sum(unnormalized), sum(posterior))## sum(unnormalized) sum(posterior)
## 1 0.8765603 1Next we’ll examine the grid approxiation posterior rounded to 2 decimal places.
#step 4
round(grid_data, 2) #2 is indicating how many decimal points to use ## pi_grid prior liklihood unnormalized posterior
## 1 0.00 0.00 0.00 0.00 0.00
## 2 0.25 3.00 0.28 0.85 0.96
## 3 0.50 0.70 0.04 0.03 0.04
## 4 0.75 0.01 0.00 0.00 0.00
## 5 1.00 0.00 0.00 0.00 0.00Finally, we’ll plot this model.
# plotting the grid approximation posterior
ggplot(grid_data, aes(x = pi_grid, y = posterior)) +
geom_point() +
geom_segment(aes(x = pi_grid, xend = pi_grid, y = 0, yend = posterior))
- Now we’re asked to repeat the above using a grid of 201 equally spaced values between 0 and 1.
#step 1
grd_data <- data.frame(p_grid = seq(from = 0, to = 1, length = 201)) #changed length to 201
#step 2
grd_data <- grd_data |> #changed names of data, prior, liklihood to keep distinct from part a
mutate(prir = dbeta(p_grid, 3,8), #prior remained the same
liklhood = dbinom(2,10, p_grid)) #successes and outcomes remained teh same
#step 3
grd_data <- grd_data |>
mutate(unnormalized = liklhood * prir,
posterir = unnormalized/sum(unnormalized))
#confirming that posterior approximation sums to 1
grd_data |> summarize(sum(unnormalized), sum(posterir))## sum(unnormalized) sum(posterir)
## 1 41.79567 1#step 4
round(grd_data, 2)## p_grid prir liklhood unnormalized posterir
## 1 0.00 0.00 0.00 0.00 0.00
## 2 0.00 0.01 0.00 0.00 0.00
## 3 0.01 0.03 0.00 0.00 0.00
## 4 0.01 0.07 0.01 0.00 0.00
## 5 0.02 0.13 0.02 0.00 0.00
## 6 0.03 0.19 0.02 0.00 0.00
## 7 0.03 0.26 0.03 0.01 0.00
## 8 0.04 0.34 0.04 0.01 0.00
## 9 0.04 0.43 0.05 0.02 0.00
## 10 0.04 0.53 0.06 0.03 0.00
## 11 0.05 0.63 0.07 0.05 0.00
## 12 0.06 0.73 0.09 0.06 0.00
## 13 0.06 0.84 0.10 0.08 0.00
## 14 0.06 0.95 0.11 0.11 0.00
## 15 0.07 1.06 0.12 0.13 0.00
## 16 0.07 1.17 0.14 0.16 0.00
## 17 0.08 1.29 0.15 0.19 0.00
## 18 0.09 1.40 0.16 0.22 0.01
## 19 0.09 1.51 0.17 0.26 0.01
## 20 0.10 1.62 0.18 0.30 0.01
## 21 0.10 1.72 0.19 0.33 0.01
## 22 0.10 1.83 0.20 0.37 0.01
## 23 0.11 1.93 0.21 0.41 0.01
## 24 0.12 2.02 0.22 0.45 0.01
## 25 0.12 2.12 0.23 0.49 0.01
## 26 0.12 2.21 0.24 0.53 0.01
## 27 0.13 2.30 0.25 0.57 0.01
## 28 0.14 2.38 0.26 0.61 0.01
## 29 0.14 2.45 0.26 0.65 0.02
## 30 0.14 2.53 0.27 0.68 0.02
## 31 0.15 2.60 0.28 0.72 0.02
## 32 0.16 2.66 0.28 0.75 0.02
## 33 0.16 2.72 0.29 0.78 0.02
## 34 0.16 2.77 0.29 0.80 0.02
## 35 0.17 2.82 0.29 0.83 0.02
## 36 0.18 2.87 0.30 0.85 0.02
## 37 0.18 2.91 0.30 0.87 0.02
## 38 0.18 2.94 0.30 0.88 0.02
## 39 0.19 2.97 0.30 0.89 0.02
## 40 0.20 3.00 0.30 0.90 0.02
## 41 0.20 3.02 0.30 0.91 0.02
## 42 0.21 3.04 0.30 0.92 0.02
## 43 0.21 3.05 0.30 0.92 0.02
## 44 0.22 3.06 0.30 0.92 0.02
## 45 0.22 3.06 0.30 0.91 0.02
## 46 0.22 3.06 0.30 0.91 0.02
## 47 0.23 3.06 0.29 0.90 0.02
## 48 0.24 3.05 0.29 0.89 0.02
## 49 0.24 3.04 0.29 0.88 0.02
## 50 0.24 3.02 0.29 0.86 0.02
## 51 0.25 3.00 0.28 0.85 0.02
## 52 0.26 2.98 0.28 0.83 0.02
## 53 0.26 2.96 0.27 0.81 0.02
## 54 0.26 2.93 0.27 0.79 0.02
## 55 0.27 2.90 0.26 0.77 0.02
## 56 0.28 2.87 0.26 0.74 0.02
## 57 0.28 2.83 0.25 0.72 0.02
## 58 0.29 2.79 0.25 0.70 0.02
## 59 0.29 2.75 0.24 0.67 0.02
## 60 0.30 2.71 0.24 0.65 0.02
## 61 0.30 2.67 0.23 0.62 0.01
## 62 0.30 2.62 0.23 0.60 0.01
## 63 0.31 2.58 0.22 0.57 0.01
## 64 0.32 2.53 0.22 0.55 0.01
## 65 0.32 2.48 0.21 0.52 0.01
## 66 0.32 2.43 0.20 0.50 0.01
## 67 0.33 2.38 0.20 0.47 0.01
## 68 0.34 2.32 0.19 0.45 0.01
## 69 0.34 2.27 0.19 0.43 0.01
## 70 0.35 2.22 0.18 0.40 0.01
## 71 0.35 2.16 0.18 0.38 0.01
## 72 0.36 2.11 0.17 0.36 0.01
## 73 0.36 2.05 0.16 0.34 0.01
## 74 0.36 2.00 0.16 0.32 0.01
## 75 0.37 1.94 0.15 0.30 0.01
## 76 0.38 1.89 0.15 0.28 0.01
## 77 0.38 1.83 0.14 0.26 0.01
## 78 0.38 1.78 0.14 0.24 0.01
## 79 0.39 1.72 0.13 0.23 0.01
## 80 0.40 1.67 0.13 0.21 0.01
## 81 0.40 1.61 0.12 0.19 0.00
## 82 0.41 1.56 0.12 0.18 0.00
## 83 0.41 1.51 0.11 0.17 0.00
## 84 0.42 1.45 0.11 0.15 0.00
## 85 0.42 1.40 0.10 0.14 0.00
## 86 0.42 1.35 0.10 0.13 0.00
## 87 0.43 1.30 0.09 0.12 0.00
## 88 0.44 1.25 0.09 0.11 0.00
## 89 0.44 1.20 0.08 0.10 0.00
## 90 0.44 1.16 0.08 0.09 0.00
## 91 0.45 1.11 0.08 0.08 0.00
## 92 0.46 1.06 0.07 0.08 0.00
## 93 0.46 1.02 0.07 0.07 0.00
## 94 0.47 0.98 0.07 0.06 0.00
## 95 0.47 0.93 0.06 0.06 0.00
## 96 0.48 0.89 0.06 0.05 0.00
## 97 0.48 0.85 0.06 0.05 0.00
## 98 0.48 0.81 0.05 0.04 0.00
## 99 0.49 0.78 0.05 0.04 0.00
## 100 0.50 0.74 0.05 0.03 0.00
## 101 0.50 0.70 0.04 0.03 0.00
## 102 0.50 0.67 0.04 0.03 0.00
## 103 0.51 0.64 0.04 0.02 0.00
## 104 0.52 0.60 0.04 0.02 0.00
## 105 0.52 0.57 0.03 0.02 0.00
## 106 0.52 0.54 0.03 0.02 0.00
## 107 0.53 0.51 0.03 0.02 0.00
## 108 0.54 0.48 0.03 0.01 0.00
## 109 0.54 0.46 0.03 0.01 0.00
## 110 0.54 0.43 0.02 0.01 0.00
## 111 0.55 0.41 0.02 0.01 0.00
## 112 0.56 0.38 0.02 0.01 0.00
## 113 0.56 0.36 0.02 0.01 0.00
## 114 0.57 0.34 0.02 0.01 0.00
## 115 0.57 0.32 0.02 0.01 0.00
## 116 0.58 0.30 0.02 0.00 0.00
## 117 0.58 0.28 0.01 0.00 0.00
## 118 0.58 0.26 0.01 0.00 0.00
## 119 0.59 0.24 0.01 0.00 0.00
## 120 0.60 0.23 0.01 0.00 0.00
## 121 0.60 0.21 0.01 0.00 0.00
## 122 0.60 0.20 0.01 0.00 0.00
## 123 0.61 0.18 0.01 0.00 0.00
## 124 0.62 0.17 0.01 0.00 0.00
## 125 0.62 0.16 0.01 0.00 0.00
## 126 0.62 0.15 0.01 0.00 0.00
## 127 0.63 0.14 0.01 0.00 0.00
## 128 0.64 0.13 0.01 0.00 0.00
## 129 0.64 0.12 0.01 0.00 0.00
## 130 0.64 0.11 0.00 0.00 0.00
## 131 0.65 0.10 0.00 0.00 0.00
## 132 0.66 0.09 0.00 0.00 0.00
## 133 0.66 0.08 0.00 0.00 0.00
## 134 0.66 0.08 0.00 0.00 0.00
## 135 0.67 0.07 0.00 0.00 0.00
## 136 0.68 0.06 0.00 0.00 0.00
## 137 0.68 0.06 0.00 0.00 0.00
## 138 0.69 0.05 0.00 0.00 0.00
## 139 0.69 0.05 0.00 0.00 0.00
## 140 0.70 0.04 0.00 0.00 0.00
## 141 0.70 0.04 0.00 0.00 0.00
## 142 0.70 0.03 0.00 0.00 0.00
## 143 0.71 0.03 0.00 0.00 0.00
## 144 0.72 0.03 0.00 0.00 0.00
## 145 0.72 0.03 0.00 0.00 0.00
## 146 0.72 0.02 0.00 0.00 0.00
## 147 0.73 0.02 0.00 0.00 0.00
## 148 0.74 0.02 0.00 0.00 0.00
## 149 0.74 0.02 0.00 0.00 0.00
## 150 0.74 0.01 0.00 0.00 0.00
## 151 0.75 0.01 0.00 0.00 0.00
## 152 0.76 0.01 0.00 0.00 0.00
## 153 0.76 0.01 0.00 0.00 0.00
## 154 0.76 0.01 0.00 0.00 0.00
## 155 0.77 0.01 0.00 0.00 0.00
## 156 0.78 0.01 0.00 0.00 0.00
## 157 0.78 0.01 0.00 0.00 0.00
## 158 0.78 0.00 0.00 0.00 0.00
## 159 0.79 0.00 0.00 0.00 0.00
## 160 0.80 0.00 0.00 0.00 0.00
## 161 0.80 0.00 0.00 0.00 0.00
## 162 0.80 0.00 0.00 0.00 0.00
## 163 0.81 0.00 0.00 0.00 0.00
## 164 0.82 0.00 0.00 0.00 0.00
## 165 0.82 0.00 0.00 0.00 0.00
## 166 0.83 0.00 0.00 0.00 0.00
## 167 0.83 0.00 0.00 0.00 0.00
## 168 0.84 0.00 0.00 0.00 0.00
## 169 0.84 0.00 0.00 0.00 0.00
## 170 0.84 0.00 0.00 0.00 0.00
## 171 0.85 0.00 0.00 0.00 0.00
## 172 0.86 0.00 0.00 0.00 0.00
## 173 0.86 0.00 0.00 0.00 0.00
## 174 0.86 0.00 0.00 0.00 0.00
## 175 0.87 0.00 0.00 0.00 0.00
## 176 0.88 0.00 0.00 0.00 0.00
## 177 0.88 0.00 0.00 0.00 0.00
## 178 0.88 0.00 0.00 0.00 0.00
## 179 0.89 0.00 0.00 0.00 0.00
## 180 0.90 0.00 0.00 0.00 0.00
## 181 0.90 0.00 0.00 0.00 0.00
## 182 0.90 0.00 0.00 0.00 0.00
## 183 0.91 0.00 0.00 0.00 0.00
## 184 0.92 0.00 0.00 0.00 0.00
## 185 0.92 0.00 0.00 0.00 0.00
## 186 0.92 0.00 0.00 0.00 0.00
## 187 0.93 0.00 0.00 0.00 0.00
## 188 0.94 0.00 0.00 0.00 0.00
## 189 0.94 0.00 0.00 0.00 0.00
## 190 0.95 0.00 0.00 0.00 0.00
## 191 0.95 0.00 0.00 0.00 0.00
## 192 0.96 0.00 0.00 0.00 0.00
## 193 0.96 0.00 0.00 0.00 0.00
## 194 0.96 0.00 0.00 0.00 0.00
## 195 0.97 0.00 0.00 0.00 0.00
## 196 0.98 0.00 0.00 0.00 0.00
## 197 0.98 0.00 0.00 0.00 0.00
## 198 0.98 0.00 0.00 0.00 0.00
## 199 0.99 0.00 0.00 0.00 0.00
## 200 1.00 0.00 0.00 0.00 0.00
## 201 1.00 0.00 0.00 0.00 0.00Next I’ll plot this approximation.
ggplot(grd_data, aes(x = p_grid, y = posterir)) +
geom_point() +
geom_segment(aes(x = p_grid, xend = p_grid, y = 0, yend = posterir))
We can see that b gives a much better estimate than a.
Exercise 6.6
We’re given a Gamma-Poison model for Y|lambda~Pois(lambda) and lambda~Gamma(20,5). N = 3 and (Y1, Y2, Y3) = (0, 1, 0).
- Use grid approximation with values E{0, 1, 2,…8} to approximate posterior model lambda.
# Step 1: Define a grid of 10 lambda values
grid_data <- data.frame(lambda_grid = seq(from = 0, to = 8, length = 8))
# Step 2: Evaluate the prior & likelihood at each lambda
grid_data <- grid_data |>
mutate(prior = dgamma(lambda_grid, 20, 5), #given prior
likelihood = dpois(0, lambda_grid) * dpois(1, lambda_grid) * dpois(0, lambda_grid))
# Step 3: Approximate the posterior
grid_data <- grid_data |>
mutate(unnormalized = likelihood * prior,
posterior = unnormalized / sum(unnormalized))
# Set the seed
set.seed(1500)
# Step 4: sample from the discretized posterior
post_sample <- sample_n(grid_data, size = 10000,
weight = posterior, replace = TRUE)
#plotting histogram of the above
ggplot(post_sample, aes(x = lambda_grid)) +
geom_histogram(aes(y = ..density..), color = "white") +
stat_function(fun = dgamma, args = list(20, 5)) +
lims(x = c(0, 8))## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.## Warning: Removed 2 rows containing missing values (geom_bar).
- Repeat the above using a grid of 201 with equally spaced values between 0 & 8.
# Step 1: Define a grid of 10 lambda values
grid_data <- data.frame(lambda_grid = seq(from = 0, to = 8, length = 201)) # went up to 10 values to see what happens after 8
# Step 2: Evaluate the prior & likelihood at each lambda
grid_data <- grid_data |>
mutate(prior = dgamma(lambda_grid, 20, 5), #given prior
likelihood = dpois(0, lambda_grid) * dpois(1, lambda_grid) * dpois(0, lambda_grid))
# Step 3: Approximate the posterior
grid_data <- grid_data |>
mutate(unnormalized = likelihood * prior,
posterior = unnormalized / sum(unnormalized))
# Set the seed
set.seed(1500)
# Step 4: sample from the discretized posterior
post_sample <- sample_n(grid_data, size = 10000,
weight = posterior, replace = TRUE)
#plotting histogram of the above
ggplot(post_sample, aes(x = lambda_grid)) +
geom_histogram(aes(y = ..density..), color = "white") +
stat_function(fun = dgamma, args = list(20, 5)) +
lims(x = c(0, 8))## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.## Warning: Removed 2 rows containing missing values (geom_bar).
## Exercise 6.7
We’re given a normal-normal model for mu with Y|mu~N(mu, 1.3^2) and mu~N(10,1.2^2). Observed data is (Y1, Y2, Y3, Y4) = (7.1, 8.9, 8.4, 8.6)
# Step 1: Define a grid of 50 mu values
grid_data <- data.frame(mu_grid = seq(from = 7, to = 10, length = 50)) #narrowed the sequence range
#changed to 50 values to demonstrate that I know as grid spacing decreases, quality of approximation decreases
# Step 2: Evaluate the prior & likelihood at each mu
grid_data <- grid_data |>
mutate(prior = dnorm(mu_grid, mean = 10, sd = 1.2),
likelihood = dnorm(7.1, mean = mu_grid, sd = 1.3)*
dnorm(8.9, mean = mu_grid, sd = 1.3)*
dnorm(8.4, mean = mu_grid, sd = 1.3)*
dnorm(8.6, mean = mu_grid, sd = 1.3))
# Step 3: Approximate the posterior
grid_data <- grid_data |>
mutate(unnormalized = likelihood * prior,
posterior = unnormalized / sum(unnormalized))
ggplot(grid_data, aes(x = mu_grid, y = posterior)) +
geom_point() +
geom_segment(aes(x = mu_grid, xend = mu_grid, y = 0, yend = posterior))
- Repeat the above using a grid of 201 equally spaced values b/w 5 and 15.
# Step 1: Define a grid of 201 mu values
grid_data <- data.frame(mu_grid = seq(from = 5, to = 15, length = 201))
# Step 2: Evaluate the prior & likelihood at each mu
grid_data <- grid_data |>
mutate(prior = dnorm(mu_grid, mean = 10, sd = 1.2),
likelihood = dnorm(7.1, mean = mu_grid, sd = 1.3)*
dnorm(8.9, mean = mu_grid, sd = 1.3)*
dnorm(8.4, mean = mu_grid, sd = 1.3)*
dnorm(8.6, mean = mu_grid, sd = 1.3))
# Step 3: Approximate the posterior
grid_data <- grid_data |>
mutate(unnormalized = likelihood * prior,
posterior = unnormalized / sum(unnormalized))
ggplot(grid_data, aes(x = mu_grid, y = posterior)) +
geom_point() +
geom_segment(aes(x = mu_grid, xend = mu_grid, y = 0, yend = posterior))
Exercise 6.8
A high dimensional Bayesian model could be useful for analyzing how well a medication works over time across different patients. The dimensions are time/age of patient and symptoms; there will be different symptoms across patients and across time.
The issue with dimensionality is that even for 2 dimensions, there could be an enormous scale of approximations to make because approximations would need to be calculated for both the x and y axis. If we had 100 values for the x-dimension and 100 for the y-dimension that would be 10,000 approximations. This may be timely and is inefficient to run as dimensions and number of values increase.
Exercise 6.9
The drawback for both grid approximations and MCMC are that they are discretized approximations. We should remember this when using these methods for continuous data.
The advantage of both MCMC and grid approximation are that they allow for approximations for posteriors that are difficult to solve by-hand.
An advantage that grid has over MCMC is that the code seems generally easier to write and can be done in R (without using Rstan and associated C++ issues).
An advantage of MCMC over grid approximation is the ability to estimate posteriors for multiple dimensions.
Exercise 6.10
For something to be a Markov chain, event i is dependent only on i-1. There is memorylessness and otherwise independence for the event for i-2, i-3, and so on.
This scenario is not an Markov chain if I assume that the first several nights in a row involved me going to a Thai restaurant. The prior frequency of Thai restaurants (on days i-2, i-3, etc.) may influence my choice on day i.
This is also not a Markov chain because there is complete independence from day to day for winning the lottery. Markov chains need i to be dependent on i-1.
This is not a Markov chain for similar reasons to a. There is memory about i-2, i-3, etc., i.e. an experienced chess player does remember the prior moves of other players. If the players are not good and can only remember the last game they played, then this could be a Markov chain.
Exercise 6.11
Use information to write RStan syntax.
bb_model <- " data { int<lower = 0, upper = 12> Y; } parameters { real<lower = 0, upper = 1> pi; } model { Y ~ binomial(20, pi); pi ~ beta(1, 1); } "
gp_model <- " data { int<lower = 0> Y[3]; } parameters { real<lower = 0> lambda; } model { Y ~ poisson(lambda); lambda ~ gamma(4, 2); } "
normal_model <- " data { vector[4] Y; } parameters { real mu; } model { Y ~ normal(mu, 1^2); mu ~ normal(0, 10^2); } "
Exercise 6.12
Write the Rstan syntax for simulating the posterior.
bb_sim <- stan(model_code = bb_model, data = list(Y = 12), chains = 4, iter = 1000*2, seed = 100) #chose 1000 to reduce processing time
gp_sim <- stan(model_code = gp_model, data = list(Y = c(3)), chains = 4, iter = 1000*2, seed = 200)
d <- list(Y = c(12.2)) nn_sim <- stan(model_code = normal_model, data = d, chains = 4, iter = 1000*2, seed = 300)